Dating the molecular clock in fungi – how close are we?

Integration of fungal evolution with the dates of plate tectonic movements, paleoecology, and the evolution of plants and animals requires a molecular clock. Imperfect though they may be, molecular clocks provide the means to convert molecular change into geological time. The relationships among clocks, phylogeography, fossils, and substitution rate variation, along with incorporation of uncertainty into clock estimates are the topics for this commentary. This commentary is timely because, for deeper divergences on the order of hundreds of millions of years, estimates of age of origin are benefiting from increasingly accurate organismal phylogenies and increasingly realistic models of molecular evolution. Taking advantage of Bayesian approaches permitting complex assumptions about node ages and molecular evolution, we used the program BEAST to apply a relaxed lognormal clock analysis to a data set comprising 50 loci for 26 taxa. In the resulting tree, branches associated with nodes calibrated by fossils showed more dramatic substitution rate variation than branches at nodes lacking calibration. As a logical extension of this result, we suspect that undetected rate variation in the uncalibrated parts of the tree is as dramatic as in the calibrated sections, underscoring the importance of fossil calibration. Fortunately, new and interesting fungal fossils are being discovered and we review some of the new discoveries that confirm the ancient origin of important taxa. To help evaluate which fossils might be useful for constraining the ages of nodes, we selected fossils thought to be early members of their clades and used ribosomal or protein-coding gene sequence substitution rates to calculate whether fossil age and expected lineage age coincide. Where ages of a fossil and the expected age of a lineage do coincide, the fossils will be particularly useful in constraining node ages in molecular clock analyses. a 2010 The British Mycological Society. Published by Elsevier Ltd. All rights reserved. 1. The negative evidence dilemma: the more earliest fossil evidence of plants on land, e.g., Ascomycota ancient the age estimate, the lower the likelihood of finding contradictory evidence Heckman et al. (2001) published age estimates for colonization of earth by fungal and plant lineages that shocked mycologists and botanists because they were about twice as old as the 0; fax: þ1 604 822 6089. ca (M. L. Berbee). ritish Mycological Society and Basidiomycota diverged 1.2 billion years ago, and mosses diverged from vascular plants 680 million years ago. Yet Heckman et al. (2001) were the first to apply data from many loci per organism to questions of dating in fungi and, under the reasonable assumption that sampling many genes would average out gene-specific selective pressures, their efforts . Published by Elsevier Ltd. All rights reserved. 2 M. L. Berbee, J. W. Taylor resulted in a data set that kept good molecular time. Their analysis also produced the huge gap between dates consistent with the fossil record and dates consistent with the molecular data, which raised these important questions for the field. If the molecular dates were right, why does the fossil record lag hundreds of millions of years behind divergence times? If the molecular dates were wrong, what misunderstanding of molecular evolution could account for the skewing of divergence time estimates? Fossils must be used to tie a molecular phylogeny to a geologic time scale, and linking fossils and molecular estimates is not without its problems. The ancient divergence times in the Heckman et al. (2001) paper illustrate a general difficulty in dating: The older the age estimate, the lower the likelihood that any conceivable geological or fossil evidence could ever contradict it. Maximum ages provided by the age of the earth (w4.6 Ba) or the development of an oxidizing atmosphere (w2 Ba) are too ancient to be of much use in limiting the age of origin of fungi. By their nature, fossils provide minimum ages for divergences and genetic lineages can be much older than even the oldest fossil representative (Fig. 1). The most useful fossils are those that represent very early members of their lineages. However, the odds of fossilization of the first members of a new group are surely small, as are the odds of finding and recognizing these fossils. Usually, the only way to assign a maximum age to an event involves setting some sort of probabilistic distribution for the likelihood of finding a fossil representing a clade, as discussed in the next section. Except for rare cases where a continuous fossil record documents a phenotypic change that is clearly associated with the emergence of a clade, e.g., the change in pollen morphology associated with the Fig. 1 – Relationship of a fossil assignable to a particular branch and the divergence (indicated by a star) leading to that branch. (A) If a fossil was formed shortly after a divergence, its geological age should agree with ages determined from divergences of DNA sequences from taxa 1 or 2 vs. taxon 3. (B) If a fossil was formed long after a divergence, its geological age should be significantly more recent than the node age estimated from nucleotide divergences between taxa 1 and 2 vs. taxon 3. emergence of eudicots or the change in skeletal morphology associated with the divergence of birds and mammals, maximum ages for divergences will always be open to dispute. Given the paucity of fungal fossils, Heckman et al. (2001) calibrated their fungal phylogeny using a date of 1600 Ma for the divergence of fungi, animals and plants. This date had been inferred from a universal phylogeny of eukaryotes that used the bird/mammal split at about 310 Ma as the calibration point (Wang et al., 1999). The assumption underlying the dates of Wang et al. (1999) and Heckman et al. (2001) was that the molecular clock ticks at the same rate throughout eukaryotes. Peterson et al. (2004) showed that it was possible to examine that assumption using multiple calibration points in a comparison of vertebrates and invertebrates to demonstrate that molecular clocks do not tick at a uniform rate. Peterson et al. (2004) also showed that Heckman et al. (2001), by picking a slowly evolving clade for calibration, arrived at very old ages for rapidly evolving lineages such as fungi and plants. Furthermore, analytical methods could not compensate for a lack of calibration points. Even with a phylogenetic method that allowed for substitution rate variation, when using the bird/ mammal divergence as a single calibration, we found surprisingly old estimated ages for all lineages within fungi (Taylor and Berbee, 2006). Demonstrating the extent to which a single calibration point can control node ages, Lücking et al. (2009) rescaled previously published fungal phylogenies to the same fossil calibration, the 400 Ma fossil fungus Paleopyrenomycites devonicus (Taylor et al., 2004, 2005; Taylor et al., 1999). Even though the original phylogenies were made using different genes and different methods, after rescaling, they showed similar node ages (Lücking et al., 2009). Consistent though those fungal node ages may be, their reliability will, however, be in doubt until rate variation within fungi is better modeled and until calibration of molecular dates is improved. Rate variation is so clearly a feature of molecular evolution that most recent approaches to molecular clocks, as reviewed by Welch and Bromham (2005) and by Lepage et al. (2007), have been directed towards allowing rates to vary across different branches in a tree. Fungal rate variation was recognized early (Berbee and Taylor, 1993). More recently, rate variation has been investigated using likelihood methods designed to accommodate rate variation, that is, to relax a strict molecular clock (Taylor and Berbee, 2006). In the next sections we further investigate variation in substitution rates in fungi by using Bayesian methods as implemented in the program BEAST ver. 1.4.8 ; (Drummond and Rambaut, 2007). After describing briefly BEAST’s Bayesian approach to rate variation and fossil calibration, we provide an example of the application of BEAST using a data set from our previous studies (Taylor and Berbee, 2006). Bayesian methods require prior knowledge (prior probability distributions) to inform the analysis and influence the resulting posterior probability (Felsenstein, 2004). As they run, Bayesian programs including BEAST launch a succession of millions of Markov chain Monte Carlo (MCMC) generations. MCMC is a trial-and-error approach where selection usually favours replacement of less likely values (given the priors) with more likely values for evolutionary parameters. However, to explore a wide range of combinations of parameter values, occasional replacements by less likely values are Dating the molecular clock in fungi 3 also permitted. With each generation, the likelihood (given the priors) of a slightly different parameter set is calculated. In running BEAST, likelihood increases with initial generations and then reaches a nearly stable plateau. The goal then is to establish posterior distributions of evolutionary parameter values in the form of samples from the generations after likelihood stabilizes. For example, a node age would be estimated as 350 Ma if that were the mean age from the sampled generations, and the confidence interval would extend from 300 to 400 Ma if 95 % of the estimated ages lay between these bounds. In this manner, at the end of the trial-and-error process, the frequencies of clades, node ages or rates for branches in the posterior distribution provide estimates of their posterior probabilities (Drummond and Rambaut, 2007). 2. Incorporating uncertainty about fossil calibrations into node age estimates: an example using the program BEAST and a Bayesian approach to a 50-gene data set To explore the consequences of a Bayesian approach to dating nodes, we used BEAST to analyze a 50-gene data set consisting of amino acids inferred from codons found in DNA sequence for 26 taxa described earlier (Rokas et al., 2005; Taylor and Berbee, 2006). Using BEAUTI (a program distributed with BEAST), we set priors for the analyses, and produced the necessary, correctly formatted XML input file for BEAST. Based on our earlier analysis (Taylor and Berbee, 2006), we used a Wagner model of amino acid substitution and a gamma site heterogeneity model with four rate categories as priors. BEAST analyses will not run if the likelihood of the combination of starting parameter values is too low. In order to run an analysis that included calibration points, we had to provide a sufficiently likely prior user tree. To create the tree, we ran BEAST without a calibration point for 500,000 generations. We then edited the XML file to include the resulting tree in Fig. 2 – Alternative options for prior probability distributions on BEAST software (Drummond and Rambaut, 2007). The x-axes re represent the instantaneous probability of the age. We used the age of divergence of bird from mammal. In this case, a 300 mil possible divergence time, and the tail of the curve is consistent place much earlier. The normal curve in (B) describes the prior animals plus fungi. The actual age is actively debated and this deviation about a prior mean of 1700 Ma. We used the uniform between 100 and 200 Ma without specifying a particular date w Newick format (a standard format, used in PHYLIP etc., see ), as a prior for further analysis. As recommended in BEAST’s documentation, we used a Yule speciation process, which specifies a constant rate of species divergence. As noted at the outset, due to the vagaries of preservation, discovery and interpretation, few fossils can be expected to accurately date a divergence. BEAST and other Bayesian programs allow the user to incorporate uncertainty about calibration ages into date estimates. Soft bounds can be applied to specify the prior probability of different ages and the final posterior age distribution is based not only on the age priors but also on the molecular substitution rates. Hard bounds, on the other hand, limit the node age estimate by controlling the range of ages considered in the analysis. For the common case where the age of a fossil provides a minimum age but the true date for a divergence may be older than the earliest representative fossil (Fig. 1), BEAST allows the user to specify an exponential distribution with a hard bound representing a minimum age from a fossil and a soft bound representing the distribution of probabilities of even older ages. After initial experimentation, we used an exponential distribution with a minimum age of 300 Ma and a standard deviation of 30 Ma for the divergence of birds from mammals (Fig. 2A). With the soft bounds provided by a normal distribution, dates older or younger than a specified date are possible, but less probable (Fig. 2B). Reflecting uncertainty about the true minimum or maximum age of the root where the plants diverged from the opisthokonts (fungi þ animals), we assigned the node age a prior normal distribution with a mean of 1700 million years and a standard deviation of 300 million years (Fig. 2B). For the fly/mosquito divergence, we used a normal prior with a mean age of 235 Ma (Peterson et al., 2004) and a standard deviation of 24 Ma. Where a continuous fossil record exists, as with the eudicots divergence, hard bounds for the minimum and even maximum age can be specified. For the divergence of eudicots (Arabidopsis) from monocots (rice), we used ages of nodes, given a fossil calibration point, as provided in present age in millions of years; the y-axes in (A), (B), exponential curve (A) to provide the prior probability for the lion-year old fossil provided a convincing lower limit for with the possibility that the divergence could have taken probability for the age of the split of the plants from the was reflected in our choice of a large (300 Ma) standard distribution (C) to constrain the age of origin of the eudicots ithin the range. 4 M. L. Berbee, J. W. Taylor a uniform prior with an upper hard bound of 200 Ma and a lower hard bound of 100 Ma (Fig. 2c). At each generation, BEAST estimated the likelihood of the data, given the priors, the tree topology and other program parameters. Lacking the ability to run BEAST in parallel, this task was computationally intensive. Using three 2.66 GHz Intel Core 2 Duo Mac Os X ver. 10.5.6 computers with 2 GB of available RAM, the 15,713,400 completed generations that we analyzed came from four independent runs of the BEAST program; one 5 million generation run that was complete after more than 27 days and three other runs that were halted before completion. We collected one tree per 100 generations. As the criterion to recognize that the number of generations was sufficient, we looked for stable likelihood plateaus from independent runs using the graphing function in Tracer (a computer program supplied with BEAST), discarding as burnin the less likely first 5001 sampled generations per run. The mean likelihoods from the four runs were reassuringly similar, ranging from ln 236,455 to 236,459. After discarding the burnin, 137,130 trees remained to be analyzed. The cumulative overall effective sample size (an estimate of the number of independent samples of the posterior distribution) was over 200, an acceptable level according to documentation in the BEAST program manual. The consensus of the 137,130 trees Fig. 3 – Amino acid substitution rates vary throughout plants, a thickness and colour throughout this phylogram. This tree was and rooted with the outgroup Plasmodium falciparum. Branch th colour indicates the rates from red for the highest substitution r rates. Only seven error bars (in black) are shown and each repre error bars were comparable in size to the five longest bars. Note t and rice, and exceptionally low rates from the node to the chick around both nodes. Narrow confidence intervals resulted from t hard bounds on the node dates, and they may be an artefact of is shown in Fig. 3. The wide error bars for node ages represented the many equally likely solutions for age estimates in the posterior distributions, even under the relaxed clock model that allowed for substitution rate changes (Fig. 3). Models of sequence evolution and soft bounds on fossil dates underestimate localized rate variation producing unreasonable node ages The appeal of the Bayesian approach lies in its ability to provide a simultaneous probabilistic estimate of values for a highly complex set of parameters including tree topology and node ages, given the priors and given the data. However, weighing against these strengths, whether or not it is possible to select reasonable priors is a serious concern in Bayesian analysis (Felsenstein, 2004). Different priors for node ages had a striking effect on our results. Given soft bounds as the priors for the age of eudicots (provided in preliminary experiments as an exponential distribution with a minimum age of 100 Ma or a normal distribution with a mean of 100 Ma) BEAST allowed the node date for the origin of eudicots to creep back in time to over 350 Ma, an unreasonable date well outside of the prior probable range (results not shown). Our large amount of sequence nimals and fungi as highlighted by changes in branch generated with the program BEAST, using four calibrations ickness is proportional to reconstructed substitution rates; ates, through purple for intermediate and blue for the lowest sents the 95 % posterior distribution for a node age. Omitted hat the exceptionally high rates from the node to Arabidopsis en and human, along with the narrow confidence intervals ension between substitution rates and fossil-based enforced the methods rather than an indication of high precision. Dating the molecular clock in fungi 5 data per taxon along with the high substitution rates in the plants overwhelmed the information content of the fossilbased prior probability distribution. Imposing a hard bound in the form of a uniform prior of 100–200 Ma, the age estimate for the node was 200 Ma, at the limit of the bound interval with unusually high substitution rates on the descendant branches from the node (Fig. 3). The narrow confidence interval around the node (Fig. 3) reflected tension between the maximum 200 Ma age from fossil calibration and the low probability of a lesser age that resulted from the model of substitution rate variation. In the absence of an imposed minimum age, the age estimate of the mammal/bird split was an impossibly recent w100 Ma (Taylor and Berbee, 2006). Forcing the mammal/bird node to a minimum age of 300 Ma (as documented through an excellent fossil record) resulted in exceptionally low substitution rates in the descendent lineages (Fig. 3). Here, the narrow confidence interval around the node reflected tension between the minimum 300 Ma age from fossil calibration and the improbability of an age this great, given the overall model of substitution rates. While these branches calibrated by hard bounds showed dramatic substitution rate variation, the uncalibrated parts of the tree may have had equally dramatic substitution rate shifts that could not be detected. In the above examples as well as in some computer simulations (Yang and Rannala, 2006), narrow confidence intervals and dates near the limits of a hard bound indicated localized trouble rather than reliable precision in dating. Yang and Rannala (2006) (Fig. 4) showed that one bad calibration point at the root of a tree, in this case the assumption that the date range for a root with a hard bound was 3.5–4.5 when it actually was 1.0, gave ages for the root that were exactly 3.5 Fig. 4 – Comparison of soft bounds and hard bounds when usin a computer simulation by Yang and Rannala (2006) shows the axis represents age and the ‘x’ axis represents the amount of s 1.0 Ma, but the calibrating fossil was ‘mistakenly’ assigned a ra dates younger than 3.5 Ma, albeit with reduced probability, so t calibration points partially compensated for the erroneous root d for the erroneously dated fossil. Disturbingly, with more data, t mating precision and underestimating error. Figure reprinted fr with permission from Oxford University Press. times too old. Even worse, at the end of analysis, the incorrect age estimates had narrow confidence intervals. Note that the analysis was precise due to the amount of sequence data, and that the most recent possible age, given the constraint, was chosen, but that it was the fossil calibration that determined the date, an inaccurate one in this case. Given a soft bound for the calibration, enough sequence data, and an adequate model of substitution, the estimate of the age of the node improved, although the apparent precision was still too high, given the true error. Variation in nucleotide substitution rates is thought to correlate with the frequency of nuclear division. For example, in vertebrates in general (Ellegren, 2007) or in mammals in particular (Fontanillas et al., 2007), substitution may be slow because the female germ line undergoes few mitotic divisions relative to generation time. If mitotic frequency correlates with overall substitution rate in fungi, where there is no separate germ line, then rate variation might be used to estimate life history parameters, such as, growth rate or the time between germination and reproduction for fungi living in substrates where direct observation is impossible. 3. Will new fossil finds help establish minimum clade ages? In this next section, we offer predictions about which recent fossil finds are old enough to increase the precision of node age estimates. The broad fungal phylogenies that we have just discussed were produced with large data sets comprising 50 genes and our Bayesian analyses showed wide variation in evolutionary rates and statistical uncertainty about node ages. g erroneous fossil divergence dates. This graph from posterior means and 95 % confidence interval where the ‘y’ equence data. In the simulation, the actual root date was nge of dates from 3.5 to 4.5 Ma. Left. Soft bounds permitted hat the combination of more data and three other good ate. Right. Hard bounds prevented the use of younger dates he confidence interval shrank to a point, greatly overestiom Yang and Rannala (2006), fig. 6, p. 220, copyright (2006), 6 M. L. Berbee, J. W. Taylor Most fungal lineages with early fossil records cannot yet be included in the broad analysis because they are represented by too few DNA sequences. For these key lineages, increasing the amount of sequence data may ultimately improve molecular dating. To evaluate the ability of specific fossils to constrain minimum ages, we simplified matters by using the available conserved ribosomal SSU gene sequences for deep divergences. From our earlier experiments with molecular clocks and under our assumption that the divergence of animals from fungi was about 1 Ba, the average rate of substitution in the nuclear small subunit ribosomal genes was about 1– 1.25 % per lineage per 100 Ma (Berbee and Taylor, 1993, 2001). Heckman et al. (2001) took the divergence of animals from fungi to be 1.6 Ba, which implied a slower rate of substitution, 0.78 %. Using our estimates of substitution in SSU per 100 Ma, a fossil is most likely to constrain the minimum age of a node if the ratio of its age to phylogenetic distance is about 100 million years for every 1.25 % substitution per lineage in the nuclear ribosomal small subunit gene. Because the rate estimates from Heckman et al. (2001) are lower, a fossil would have to be older relative to a node to constrain an age. By any estimate, if the substitution percentage between sister taxa is much higher than 2 % per lineage per 100 Ma when calibrated by the fossil, the lineage is likely much older that the fossil (although it could also be evolving unusually quickly). We (i) arbitrarily selected two species with roughly average branch lengths that diverged from the node constrained by a specific fossil, (ii) estimated the percent substitution in SSU rDNA between the two species (iii) divided the pairwise distance by two for an estimate of the percent substitution per lineage and (iv) multiplied the per lineage percent substitution by the estimated SSU substitution rate per (100 Ma/1.25 % or 100 Ma/2 %) to estimate the time since divergence. We then compared the age of the fossil with the estimated time since divergence. For the estimates of pairwise distance, we used a Kimura two-parameter (K2P) model if the percent substitution was less than 5 %, and a general time reversible (GTR) model of substitution, with rate variation described by a gamma distribution using a likelihood estimate for a gamma parameter, if the percent substitution was over 5 %. For recent nodes where change in the conserved SSU has been too slow to permit useful estimates, the more rapid substitutions at neutral sites in protein-coding genes served to test if a fossil was likely to provide a meaningful constraint. Neutral rates can be estimated from any protein-coding gene DNA sequence where a substitution does not change an amino acid because these synonymous nucleotide substitutions are thought to elude the forces of natural selection. Although factors such as codon bias can affect rates of synonymous substitution, they can be ignored for our rough estimates of fossil utility. Neutral rates of protein gene evolution in fungi, plants, animals and bacteria are typically between 1 10 8 and 10 9 substitutions per site per year (Kasuga et al., 2002) and so we divided the per lineage percent substitution by 1 10 9 to estimate node age in years for comparisons with fossil age. Fossil Glomeromycota This group of fungi, generally accepted as the phylogenetic sister clade to the Ascomycota and Basidiomycota, consists of arbuscular mycorrhizal symbionts. The Glomeromycota have left excellent and ancient fossils. Redecker et al. (2000) described fossil spores from shallow marine sediments from the Ordovician, 455–460 Ma in age, that closely resemble Glomeromycota spores. Arbuscules, finely-branched Glomeromycota fungal hyphae within plant cells, were clearly preserved in cells of stems of Aglaophyton, a primitive land plant from 400 Ma Rhynie chert (Remy et al., 1994) and from roots from the Triassic (250–199 Ma) (Stubblefield et al., 1987). The pairwise difference between Glomus eburneum SSU sequence, GenBank accession number (GB) AM713431 and either a basidiomycete (Agaricus bisporus GB L36658) or an ascomycete (Peziza gerardii GB DQ646543) is about 19 %, or about 9.5 % per lineage (GTR model). These fossils are tremendously important and they provide some of the best and most readily interpreted evidence for the presence of their phylum. At a nucleotide substitution rate of 1.25 %/100 Ma, the oldest fossil at 460 Ma would be far younger than the divergence at 760 Ma. It is not until the substitution rate approaches 2 %/100 Ma that the divergence, at 475 Ma, would approach the fossil age. It seems likely that the Glomeromycota stem lineage arose earlier than suggested by the fossils and that it predated land plants. Whatever characters or features of habitat predisposed the ancient Glomeromycota towards symbiosis remain a mystery.